require(estimatr)

#### conjoint data ####
conjoint.data <- read.csv("replication_data_conjoint.csv")

# number of respondents
max(conjoint.data$respondent.id)
# number of observations
nrow(conjoint.data)

# reorder factor levels
conjoint.data$party <- factor(conjoint.data$party, 
                              levels = c("Independent", "LDP", "CDP", 
                                         "Komeito", "JIP", "JCP"))
conjoint.data$cand.gender <- factor(conjoint.data$cand.gender, 
                                    levels = c("Male", "Female"))
conjoint.data$cand.age <- factor(conjoint.data$cand.age, 
                                 levels = c("42", "52", "59", "67"))
conjoint.data$cand.edu <- factor(conjoint.data$cand.edu, 
                                 levels = c("High school", "Vocational college", 
                                            "Private university", "Local public university", 
                                            "The University of Tokyo", "Graduate school"))
conjoint.data$occupation <- factor(conjoint.data$occupation, 
                                   levels = c("Buniness employee", "Business executive", 
                                              "Government employee", "Celebrity", 
                                              "Secretary", "Local politician"))
conjoint.data$hometown <- factor(conjoint.data$hometown, 
                                 levels = c("Inside", "Outside"))
conjoint.data$experience <- factor(conjoint.data$experience, 
                                   levels = c("No experience", "6 years in office", 
                                              "12+ years in office"))
conjoint.data$legacy <- factor(conjoint.data$legacy, 
                               levels = c("None", "Local politician", 
                                          "Diet member", "Cabinet minister"))

# respondents' gender category
conjoint.data$gender.cat <- ifelse(conjoint.data$gender == 3, NA, 
                                   conjoint.data$gender)
# respondents' age category
conjoint.data$age.cat <- ifelse(conjoint.data$age < 40, 1, 
                                ifelse(conjoint.data$age < 60, 2, 3))
# respondents' education category
conjoint.data$education <- ifelse(conjoint.data$education < 3, 1, 
                                  ifelse(conjoint.data$education < 6, 2, 3))
# respondents' knowledge category
conjoint.data$knowledge.cat <- ifelse(conjoint.data$knowledge > 3, 1, 
                                      ifelse(conjoint.data$knowledge == 3, 2, 3))
# respondents' LDP support
conjoint.data$LDP.support <- (conjoint.data$partisanship == 1) * 1

#### AMCE estimation ####
ACME.results <- lm_robust(rating ~ legacy + cand.gender + cand.age + 
                            cand.edu + occupation + hometown + 
                            experience + party, 
                          data = conjoint.data, clusters = respondent.id)
round(summary(ACME.results)$coefficients, 2)

cairo_pdf("Figure_5.pdf", 
          width = 4.5, height = 1.8, pointsize = 9, family = "Helvetica")
par(mar = c(3.5, 0, 0.5, 0.5), lwd = 0.5)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-0.8, 0), ylim = c(0, 4.5), 
     xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = 0, col = "gray50")
abline(v = seq(-0.5, -0.1, 0.1), col = "gray50", lty = 3)
points(ACME.results$coefficients[2:4], 3:1, pch = 19)
segments(ACME.results$conf.low[2:4], 3:1, 
         ACME.results$conf.high[2:4], 3:1)
text(-0.8, 4, labels = "His or her parent was ...", pos = 4)
text(rep(-0.74, 3), 3:1, 
     labels = c("a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
axis(1, at = seq(-0.5, 0, 0.1), lwd = 0.5)
mtext("Average marginal component effect", side = 1, at = -0.25, line = 2.5)
dev.off()

attribute.labels <- c("Legacy status", "Gender", "Age", "Education", 
                      "Occupation", "Hometown", "Political experience", "Party")

cairo_pdf("Figure_A10.pdf", 
          width = 5, height = 6, pointsize = 8, family = "Helvetica")
par(mar = c(3.5, 0, 0.5, 0.5), lwd = 0.5)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-1.8, 0.5), ylim = c(0, 40), 
     xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = 0, col = "gray50")
abline(v = c(seq(-1, -0.25, 0.25), 0.25, 0.5), lty = 3, col = "gray50")
x <- 2
y <- 40
for (i in 1:8) {
  text(-1.8, y, attribute.labels[i], pos = 4)
  text(-1.8, y - 1, 
       paste0("(base: ", tolower(ACME.results$xlevels[[i]][1]), ")"), pos = 4)
  for (j in 1:(length(ACME.results$xlevels[[i]]) - 1)) {
    segments(-1.05, y - 1 - j, 0.55, y - 1 - j, col = "gray80")
    text(-1.07, y - 1 - j, ACME.results$xlevels[[i]][j + 1], pos = 2)
    points(ACME.results$coefficients[x], y - 1 - j, pch = 19)
    segments(ACME.results$conf.low[x], y - 1 - j, 
             ACME.results$conf.high[x], y - 1 - j)
    x <- x + 1
  }
  y <- y - 2 - j
}
axis(1, at = seq(-1, 0.5, 0.25), lwd = 0.5)
mtext("Average marginal component effect", side = 1, at = -0.25, line = 2.5)
dev.off()

#### heterogeneous treatment effects (Appendix E.3) ####
# function to estimate marginal means (MMs)
MM <- function(data, attribute, outcome, id) {
  lm.data <- data.frame(id = data[, id], 
                        A = data[, attribute], 
                        Y = data[, outcome])
  l <- nlevels(lm.data$A)
  result.matrix <- matrix(NA, l, 3)
  for (i in 1:l) {
    level.label <- levels(lm.data$A)[i]
    result.matrix[i, ] <- unlist(lm_robust(Y ~ 1, data = lm.data, 
                                           subset = A == level.label, 
                                           clusters = id)[c(1, 6, 7)])
  }
  result.matrix - mean(lm.data$Y)
}

## respondents' characteristics
# male respondents
MM.male <- MM(subset(conjoint.data, gender == 1), 
              "legacy", "rating", "respondent.id")
# female respondents
MM.female <- MM(subset(conjoint.data, gender == 2), 
                "legacy", "rating", "respondent.id")

# young respondents
MM.young.age <- MM(subset(conjoint.data, age.cat == 1), 
                   "legacy", "rating", "respondent.id")
# middle-aged respondents
MM.middle.age <- MM(subset(conjoint.data, age.cat == 2), 
                    "legacy", "rating", "respondent.id")
# older respondents
MM.old.age <- MM(subset(conjoint.data, age.cat == 3), 
                 "legacy", "rating", "respondent.id")

# respondents with low education
MM.low.edu <- MM(subset(conjoint.data, education == 1), 
                 "legacy", "rating", "respondent.id")
# respondents with middle education
MM.middle.edu <- MM(subset(conjoint.data, education == 2), 
                    "legacy", "rating", "respondent.id")
# respondents with higher education
MM.higher.edu <- MM(subset(conjoint.data, education == 3), 
                  "legacy", "rating", "respondent.id")

# respondents with low knowledge
MM.low.knowledge <- MM(subset(conjoint.data, knowledge.cat == 1), 
                       "legacy", "rating", "respondent.id")
# respondents with middle knowledge
MM.middle.knowledge <- MM(subset(conjoint.data, knowledge.cat == 2), 
                          "legacy", "rating", "respondent.id")
# respondents with high knowledge
MM.high.knowledge <- MM(subset(conjoint.data, knowledge.cat == 3), 
                        "legacy", "rating", "respondent.id")

# LDP supporters
MM.LDP <- MM(subset(conjoint.data, LDP.support == 1), 
             "legacy", "rating", "respondent.id")
# non LDP supporters
MM.nonLDP <- MM(subset(conjoint.data, LDP.support == 0), 
                "legacy", "rating", "respondent.id")

cairo_pdf("Figure_A11.pdf", 
          width = 5, height = 8, pointsize = 9, family = "Helvetica")
layout(matrix(1:5, 5, 1))
par(mar = c(3.5, 0, 3, 0.5), lwd = 0.5)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-1, 0.3), ylim = c(0, 5.5), 
     main = "Gender", xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = seq(-0.7, 0.3, 0.1), col = "gray50", lty = 3)
segments(MM.male[, 2], 4:1 + 0.1, 
         MM.male[, 3], 4:1 + 0.1)
segments(MM.female[, 2], 4:1 - 0.1, 
         MM.female[, 3], 4:1 - 0.1)
points(MM.male[, 1], 4:1 + 0.1, pch = 19)
points(MM.female[, 1], 4:1 - 0.1, pch = 21, bg = "white")
text(-1, 5, labels = "His or her parent was ...", pos = 4)
text(rep(-0.95, 4), 4:1, 
     labels = c("not a politician", 
                "a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
text(MM.female[1, 3], 4.2, "Male", pos = 4)
text(MM.female[1, 3], 3.8, "Female", pos = 4)
axis(1, at = seq(-0.7, 0.3, 0.1), 
     labels = c("-0.7", "-0.6", "-0.5", "-0.4", "-0.3", 
                "-0.2", "-0.1", "0.0", "0.1", "0.2", "0.3"), lwd = 0.5)
mtext("Centered marginal means", side = 1, 
      at = -0.25, line = 2.5, cex = 0.75)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-1, 0.3), ylim = c(0, 5.5), 
     main = "Age", xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = seq(-0.7, 0.3, 0.1), col = "gray50", lty = 3)
segments(MM.young.age[, 2], 4:1 + 0.2, 
         MM.young.age[, 3], 4:1 + 0.2)
segments(MM.middle.age[, 2], 4:1, 
         MM.middle.age[, 3], 4:1)
segments(MM.old.age[, 2], 4:1 - 0.2, 
         MM.old.age[, 3], 4:1 - 0.2)
points(MM.young.age[, 1], 4:1 + 0.2, pch = 19)
points(MM.middle.age[, 1], 4:1, pch = 21, bg = "white")
points(MM.old.age[, 1], 4:1 - 0.2, pch = 4)
text(-1, 5, labels = "His or her parent was ...", pos = 4)
text(rep(-0.95, 4), 4:1, 
     labels = c("not a politician", 
                "a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
text(MM.middle.age[1, 3], 4.4, "Young", pos = 4)
text(MM.middle.age[1, 3], 4, "Middle", pos = 4)
text(MM.middle.age[1, 3], 3.6, "Older", pos = 4)
axis(1, at = seq(-0.7, 0.3, 0.1), 
     labels = c("-0.7", "-0.6", "-0.5", "-0.4", "-0.3", 
                "-0.2", "-0.1", "0.0", "0.1", "0.2", "0.3"), lwd = 0.5)
mtext("Centered marginal means", side = 1, 
      at = -0.25, line = 2.5, cex = 0.75)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-1, 0.3), ylim = c(0, 5.5), 
     main = "Education", xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = seq(-0.7, 0.3, 0.1), col = "gray50", lty = 3)
segments(MM.low.edu[, 2], 4:1 + 0.2, 
         MM.low.edu[, 3], 4:1 + 0.2)
segments(MM.middle.edu[, 2], 4:1, 
         MM.middle.edu[, 3], 4:1)
segments(MM.higher.edu[, 2], 4:1 - 0.2, 
         MM.higher.edu[, 3], 4:1 - 0.2)
points(MM.low.edu[, 1], 4:1 + 0.2, pch = 19)
points(MM.middle.edu[, 1], 4:1, pch = 21, bg = "white")
points(MM.higher.edu[, 1], 4:1 - 0.2, pch = 4)
text(-1, 5, labels = "His or her parent was ...", pos = 4)
text(rep(-0.95, 4), 4:1, 
     labels = c("not a politician", 
                "a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
text(MM.middle.edu[1, 3], 4.4, "Low", pos = 4)
text(MM.middle.edu[1, 3], 4, "Middle", pos = 4)
text(MM.middle.edu[1, 3], 3.6, "Higher", pos = 4)
axis(1, at = seq(-0.7, 0.3, 0.1), 
     labels = c("-0.7", "-0.6", "-0.5", "-0.4", "-0.3", 
                "-0.2", "-0.1", "0.0", "0.1", "0.2", "0.3"), lwd = 0.5)
mtext("Centered marginal means", side = 1, 
      at = -0.25, line = 2.5, cex = 0.75)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-1, 0.3), ylim = c(0, 5.5), 
     main = "Political knowledge", xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = seq(-0.7, 0.3, 0.1), col = "gray50", lty = 3)
segments(MM.low.knowledge[, 2], 4:1 + 0.2, 
         MM.low.knowledge[, 3], 4:1 + 0.2)
segments(MM.middle.knowledge[, 2], 4:1, 
         MM.middle.knowledge[, 3], 4:1)
segments(MM.high.knowledge[, 2], 4:1 - 0.2, 
         MM.high.knowledge[, 3], 4:1 - 0.2)
points(MM.low.knowledge[, 1], 4:1 + 0.2, pch = 19)
points(MM.middle.knowledge[, 1], 4:1, pch = 21, bg = "white")
points(MM.high.knowledge[, 1], 4:1 - 0.2, pch = 4)
text(-1, 5, labels = "His or her parent was ...", pos = 4)
text(rep(-0.95, 4), 4:1, 
     labels = c("not a politician", 
                "a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
text(MM.high.knowledge[1, 3], 4.4, "Low", pos = 4)
text(MM.high.knowledge[1, 3], 4, "Middle", pos = 4)
text(MM.high.knowledge[1, 3], 3.6, "High", pos = 4)
axis(1, at = seq(-0.7, 0.3, 0.1), 
     labels = c("-0.7", "-0.6", "-0.5", "-0.4", "-0.3", 
                "-0.2", "-0.1", "0.0", "0.1", "0.2", "0.3"), lwd = 0.5)
mtext("Centered marginal means", side = 1, 
      at = -0.25, line = 2.5, cex = 0.75)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-1, 0.3), ylim = c(0, 5.5), 
     main = "Partisanship", xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = seq(-0.7, 0.3, 0.1), col = "gray50", lty = 3)
segments(MM.LDP[, 2], 4:1 + 0.1, 
         MM.LDP[, 3], 4:1 + 0.1)
segments(MM.nonLDP[, 2], 4:1 - 0.1, 
         MM.nonLDP[, 3], 4:1 - 0.1)
points(MM.LDP[, 1], 4:1 + 0.1, pch = 19)
points(MM.nonLDP[, 1], 4:1 - 0.1, pch = 21, bg = "white")
text(-1, 5, labels = "His or her parent was ...", pos = 4)
text(rep(-0.95, 4), 4:1, 
     labels = c("not a politician", 
                "a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
text(MM.LDP[1, 3], 4.2, "LDP supporters", pos = 4)
text(MM.LDP[1, 3], 3.8, "not LDP supporters", pos = 4)
axis(1, at = seq(-0.7, 0.3, 0.1), 
     labels = c("-0.7", "-0.6", "-0.5", "-0.4", "-0.3", 
                "-0.2", "-0.1", "0.0", "0.1", "0.2", "0.3"), lwd = 0.5)
mtext("Centered marginal means", side = 1, 
      at = -0.25, line = 2.5, cex = 0.75)
dev.off()

## election types
# HoR election
MM.HOR <- MM(subset(conjoint.data, HOC == 0), 
             "legacy", "rating", "respondent.id")
# HoC election
MM.HOC <- MM(subset(conjoint.data, HOC == 1), 
             "legacy", "rating", "respondent.id")

cairo_pdf("Figure_A12.pdf", 
          width = 5, height = 1.6, pointsize = 7, family = "Helvetica")
par(mar = c(3.5, 0, 0.5, 0.5), lwd = 0.5)
plot(NULL, NULL, type = "n", bty = "n", xlim = c(-0.75, 0.2), ylim = c(0, 5.5), 
     main = "", xlab = "", ylab = "", xaxt = "n", yaxt = "n")
abline(v = seq(-0.5, 0.2, 0.1), col = "gray50", lty = 3)
segments(MM.HOR[, 2], 4:1 + 0.1, 
         MM.HOR[, 3], 4:1 + 0.1)
segments(MM.HOC[, 2], 4:1 - 0.1, 
         MM.HOC[, 3], 4:1 - 0.1)
points(MM.HOR[, 1], 4:1 + 0.1, pch = 19)
points(MM.HOC[, 1], 4:1 - 0.1, pch = 21, bg = "white")
text(-0.75, 5, labels = "His or her parent was ...", pos = 4)
text(rep(-0.71, 4), 4:1, 
     labels = c("not a politician", 
                "a local politician", 
                "a Diet member", 
                "a cabinet minister"), 
     pos = 4)
text(MM.HOC[1, 3], 4.2, "HoR", pos = 4)
text(MM.HOC[1, 3], 3.8, "HoC", pos = 4)
axis(1, at = seq(-0.5, 0.2, 0.1), 
     labels = c("-0.5", "-0.4", "-0.3", "-0.2", 
                "-0.1", "0.0", "0.1", "0.2"), lwd = 0.5)
mtext("Centered marginal means", side = 1, 
      at = -0.25, line = 2.5)
dev.off()

## statistical tests for heterogeneous treatment effects
# function to conduct the joint significance F-test
interaction.F.test <- function(attr.name, covariate, data) {
  lm.model.1 <- lm(paste0("rating ~ legacy + cand.gender + cand.age + 
                             cand.edu + occupation + hometown + 
                             experience + party + ", covariate), data = data)
  lm.model.2 <- lm(paste0("rating ~ legacy + cand.gender + cand.age + 
                             cand.edu + occupation + hometown + 
                             experience + party + ", covariate, 
                          "+", attr.name, " : ", covariate), data = data)
  anova(lm.model.1, lm.model.2)$"Pr(>F)"[2]
}

# respondents' gender
round(interaction.F.test("legacy", "gender.cat", conjoint.data), 3)
# respondents' age
round(interaction.F.test("legacy", "age.cat", conjoint.data), 3)
# respondents' education
round(interaction.F.test("legacy", "education", conjoint.data), 3)
# respondents' knowledge
round(interaction.F.test("legacy", "knowledge", conjoint.data), 3)
# respondents' partisanship
round(interaction.F.test("legacy", "LDP.support", conjoint.data), 3)
# election type
round(interaction.F.test("legacy", "HOC", conjoint.data), 3)